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Despite having the simplest atomic structure, bulk FeSe has an observed electronic structure 
with the largest deviation from the band theory predictions among all Fe-based superconductors 
and exhibits a low temperature nematic electronic state without intervening magnetic order. We 
show that the Fe-Fe interatomic Coulomb repulsion V offers a natural explanation for the puzzling 
electron correlation effects in FeSe superconductors. It produces a strongly renormalized low-energy 
band structure where the van Hove singularity sits remarkably close to Fermi level in the high- 
temperature electron liquid phase as observed experimentally. This proximity enables the quantum 
fluctuations in V to induce a rotational symmetry breaking electronic bond order in the d- wave 
channel. We argue that this emergent low-temperature d- wave bond nematic state, different from 
the commonly discussed ferro-orbital order and spin-nematicity, has been observed recently by sev¬ 
eral angle resolved photoemission experiments detecting the lifting of the band degeneracies at high 
symmetry points in the Brillouin zone. We present a symmetry analysis of the space group and 
identify the hidden antiunitary T-symmetry that protects the band degeneracy and the electronic 
order/interaction that can break the symmetry and lift the degeneracy. We show that the d- wave 
nematic bond order, together with the spin-orbit coupling, provide a unique explanation of the tem¬ 
perature dependence, momentum space anisotropy, and domain effects observed experimentally. We 
discuss the implications of our findings on the structural transition, the absence of magnetic order, 
and the intricate competition between nematicity and superconductivity in FeSe superconductors. 

PACS numbers: 


I. INTRODUCTION 

The electron nematic phase with purely rotational 
symmetry breaking is arguably the most unconventional 
and poorly understood phase in Fe-based superconduc¬ 
tors. In the Fe-pnictides, the direct observation of ne¬ 
matic electronic structure has been difficult since the or¬ 
thorhombic lattice distortion is immediately followed by 
the collinear spin density wave (SDW) order that breaks, 
in addition to the spatial-orbital rotational symmetry, 
lattice translation, spin-rotation, and time-reversal sym¬ 
metries. This leaves the origin of nematicity highly de¬ 
bated QJ between the spin-nematic @,Hj and ferro-orbital 
order scenarios. In contrast, the Fe-chalcogenide 
FeSe undergoes the tetragonal to orthorhombic structural 
transition at T s = 87 K and the SC transition at T c = 9 K 
without any trace of magnetic order 0. The lat¬ 
ter enabled direct observations of a rotational symmetry 
breaking electronic structure at low temperatures by sev¬ 
eral angle resolved photoemission spectroscopy (ARPES) 
experiments recently (nHH , unveiling that the SC tran¬ 
sition in bulk FeSe takes place from a highly unconven¬ 
tional nematic electronic state. Understanding the mi¬ 
croscopic origin of this nematic order, its relation to the 
structural transition and magnetism is the focus of this 
work. 

In a nutshell, ARPES detects the splitting of sym¬ 
metry protected degeneracies between the d xz and d yz 
orbitals in the band dispersions at the high symmetry 
points M(7 t, 0,0) and the T(0, 0,0)/Z(0,0,7r) in the orig¬ 
inal Brillouin zone (BZ). More importantly, the corre¬ 
sponding degeneracy splitting energy Am and Ar /z ap¬ 


pear anisotropic in momentum space with different tem¬ 
perature dependence. The nematic transition is deter¬ 
mined by the strongly T-dependent Am whose onset co¬ 
incides with [I 2 I , flil o r is about 20K above the struc¬ 
tural transition T s [llL 0. Am rises with decreasing 
T and reaches Am — 62meV at 22K, closely resembling 
the T-dependence of an energy scale associated with a 
symmetry breaking order parameter. This is consistent 
with NMR 0,0|, optics, and transport measurements 
mm detecting changes in the electronic state near or 
above T s . The splitting at the BZ center, Ar — 30meV 
at 22K, on the other hand, does not break rotational 
symmetry above T s . It is nearly T -independent up to 
150K and insensitive to T s and the onset of Am 0. It 
was thus conjectured 0] that the nematic order in FeSe 
is not due to the commonly discussed ferro-orbital order 
00J23, but rather driven by a d-wave nematic bond 
order m OdNB — S/c( cos kx cos k y )[n xz (k) -\- n yz (*)]■ 
We show in this work that the nearest neighbor Fe-Fe 
interatomic Coulomb repulsion V can be the microscopic 
origin for the emergent nematic order, the absence of 
magnetism, and at a more fundamental level, the unusu¬ 
ally large band renormalization in bulk FeSe. The ob¬ 
served electronic structure of FeSe shows the largest de¬ 
viation from the local density approximation (LDA) band 
dispersions among all Fe-based superconductors 0,0 . 
A remarkable difference from the Fe-pnictides is that the 
renormalization is the strongest at low energies near the 
Fermi level (Ep) as shown in Fig. 1, where the LDA 
bands are compared with ARPES measurements. Of cru¬ 
cial importance is the symmetry protected van Hove sin¬ 
gularity (vHS) at M point created by the saddle point in 
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FIG. 1: (a) LDA band structure of FeSe.(b) Band dispersions 
observed by ARPES in the symmetric phase at 120K. (c) 
Band degeneracy at high symmetry points (F, M ) among d xz>e 
(blue), d yz , e (red), d xz>0 (cyan), and d yZ}0 (green), (d) ARPES 
results in the nematic state at 22K. Red dashed lines near T 
are data taken at 120K, showing Ar is nearly T-independent. 
Solid and dashed blue lines near M correspond to the two 
domains. 


the band dispersions of the d xz / yz orbital. Being more 
than 250meV below Ep in LDA, it moves to a mere 
E v h — 25meV below the Fermi level in ARPES in the 
high temperature electron liquid phase at T > T s . We 
show that this correlation-induced proximity of the vHS 
to Ep , concomitant with the large mass enhancement 
and the contraction of the FS pockets, sets the stage 
for the electronic nematic transition near T s that would 
ultimately gap out the vHS and induce the structural 
transition. 

What is the microscopic interaction that would drive 
the vHS so close to Ep in the high-temperature elec¬ 
tron liquid phase? The more than 10 times reduction 
in the distance from the hole band top at T to the vHS 
at M (see Fig.l) cannot come from crystal field correc¬ 
tions induced by local interactions [2f| [2g since it comes 
from the same atomic orbital. Recent LDA+DMFT (dy¬ 
namical meanfield theory) calculations [13, HH show that 
the intra-atomic Hubbard U and Hund’s rule coupling J 
would only produce a bandwidth reduction (albeit with 
some orbital dependence) and a mass enhancement by a 
factor of 3 ~ 4 and leave the FS pockets much larger than 
those observed by ARPES and quantum oscillations [l3j] . 
Thus, intra-atomic correlations alone cannot account for 
the electronic structure near Ep [2l| . 


We find that the nearest neighbor Coulomb V gener¬ 
ates directly hopping corrections to the band dispersion: 
it pushes the vHS at M up toward Ep and pulls down 
the top of the hole band at P, resulting in a low en¬ 
ergy dispersion consistent with experiments. Both the 
FS pockets and the E v h in ARPES can be produced 
by a V ~ 0.73eV in the Hartree-Fock (HF) theory with 
the bare LDA bandwidth W — 4.2eV. When the vHS is 
driven to Ep (i.e. .E v h = 0), a nematic instability of the 
Pomeranchuk-type occurs in the symmetry breaking va¬ 
lence bond channels. The leading instability corresponds 
precisely to the d-wave nematic bond order. The prox¬ 
imity to this nematic instability in bulk FeSe allows us 
to carry out a weak-coupling analysis of the extended 
t-U-V Hubbard model and find good agreement with ex¬ 
periments. Moreover, we find that the V-renormalized 
low-energy band structure with the reduced FS pockets 
promotes the d -wave nematicity while suppressing the 
collinear SDW, which is a possible explanation for the 
absence of magnetic order in bulk FeSe. 


Since the symmetry content of the band degeneracies 
and its relation to the nematic order in Fe-based su¬ 
perconductors have not been understood, we begin in 
Section II with a systematic symmetry analysis. We 
show that the band degeneracy at the high symmetry 
point originates from the existence of two-dimensional 
irreducible representations of the space gr oup involving 
rotation, reflection, and glide symmetries [23 . We then 
show that there exists two sets of “hidden” antiunitary 
T-symmetries (borrowing the T from an analogy to the 
time-reversal symmetry) that protect the band degen¬ 
eracies at the F and M points respectively. All possible 
degeneracy lifting interactions are then studied accord¬ 
ing to their symmetry and symmetry-breaking properties 
in connection to the experimental findings on the tem¬ 
perature dependence, momentum space anisotropy, and 
domain effects. The analysis shows that the only inter¬ 
actions consistent with the experimental findings are the 
atomic spin-orbit coupling governing the lifting of the 
degeneracy at F without breaking the four-fold rotation 
symmetry, and the rotational symmetry breaking d -wave 
nematic bond order that splits the degeneracy at M. In 
Section III, we develop the microscopic theory for the 
nematic state in bulk FeSe based on the t-U-V Hubbard 
model and show that the quantum fluctuations in the in¬ 
tersite correlation V leads to the important band renor¬ 
malization and in particular to the dynamical V-driven 
proximity of the vHS to the Fermi level. The nematic in¬ 
stability is studied in detail and the obtained low energy 
band structure and the Fermi surfaces are compared to 
recent experimental results. In Section IV, we provide 
a summary and discuss the implications of these find¬ 
ings on FeSe films, the effects of electron doping, and the 
interplay between nematicity and superconductivity, as 
well as propose experimental tests for the theory. 
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II. BAND DEGENERACY, SYMMETRY 
CONTENT, AND EFFECTIVE INTERACTIONS 

We begin with a discussion on the symmetry-protected 
band degeneracies at the T and M points illustrated 
in Fig. 1(c). These symmetry properties dictate a rich 
and interesting set of possible degeneracy lifting interac¬ 
tions with different implications on the momentum space 
anisotropy and the domain effects. We will show that 
only the d- wave nematic bond order is naturally consis¬ 
tent with the experimental findings for the nematic tran¬ 
sition. 

The LDA electronic structure (Fig.la) can be de¬ 
scribed by the tight-binding (TB) model down-folded to 
the Fe 3d-manifold mm, 

Ht = ^ 2 , tjfdlotcrdjpo + (1) 

crijotfi aioL 

where e a is the on-site (crystal field) energy of an electron 
in orbital a = (1, 2,3,4, 5) = ( xz , yz, x 2 — y 2 , xy, z 2 ) and 
is the electron hopping integral between sites (i,j) 
and orbitals (a, (3). The FeSe lattice structure contains 
two Fe atoms per unit cell labeled by l = A and B. As 
a result, the TB Hamiltonian in momentum space is, 


The less discussed point group symmetry E^ is a mirror 
reflection about the diagonal xy-z plane and operates in 
the spatial and orbital space according to, 

E d : x -A y,y -A x, z -A z, 

rj£ _i. rj^ _V r# A 1 _i. fit ( K 'l 

u xz 7 u yzT u yz 7 u xzi u xy 7 u xyi \ u J 

r\i fit r\t v yZ 

a x 2_ y 2 -a -a x 2_ y2 ,a z 2 -a a z2 . 

We will show that E^ plays as important a role as S 4 
in the origin of the band degeneracy. Finally, the glide 
symmetry [32j of the space group G s = E z • T x ^ y is asso¬ 
ciated with the translation by one lattice spacing along 
the x or y direction (T x , y ) followed by a mirror reflection 
about the x-y plane. Under the G s = E z • T x operation, 
A -(A B and 

G s : x —> x + a,y -A y, z -A —2: 

d xz -A -d xz ,d yz -A ~d yz ,d xy -A d xy , (6) 
r\A , j B jA , j B 

a X 2 —y 2 -A a x 2 _ y 2 ,a z 2 ^ a z 2 . 

Since R is a symmetry of H t (k ), [H t (k),R] = 0. Thus, 
R\n,k) = | n',k') is a simultaneous eigenstate of | n,k) 
at the same energy E% as given in Eq.Q. Note that in 
general k! = Rk ^ k. In particular, using Eqs. (1115161) . it 
is simple to show that in momentum space: 


H t = Y,ti(k)Ht(k)Mk), (2) 

ka 

where the basis vector ip a (,k) = [d£(k) , d^ (fc)] T with 
d e a = {d\ a ,d 2 a , ■ ■ ■ , df CT ), and H t (k) is a 10 x 10 matrix 
whose eigenvalues govern the dispersion of the 10 LDA 
bands, 

H t (k)\n,k) = E k \n,k), n = 1... 10. (3) 

The explicit form of H t {k ) and its parameters were de¬ 
rived in Refs.[23l.l3lj for FeSe. Note that to avoid confu¬ 
sion, we use the k values confined to the reduced two-Fe 
BZ shown in Fig.2(b) to label the momentum eigenstates 
such that the four M points are located at (±7r, 0) and 
(0, ±7r) of the original one-Fe BZ. 


A. Symmetry protected band degeneracies 


*84 ■ k x ^ k y , ky y k x 

E^ : k x -A k y , k y -A k x (7) 

Gs ■ k x -A k x ,k y -A k y 


However, at certain high symmetry points in the BZ (e.g. 
r and M), it is possible to have k' equal or equivalent 
to k under the reciprocal lattice vector, resulting in sym¬ 
metry related band degeneracies. To elucidate this, we 
transform the basis vector in Eq. © into the eigenbasis 
of the glide symmetry operator G s via a unitary rota¬ 
tion A’ a (k) = Uip cr (k) = U[d^(k),d^(k)) T [32|, where 



Ui U 2 
tA -U 2 
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The band degeneracies at T and M are related to three 
important space group symmetries R = S 4 , E<j, G s of the 
atomic and electronic structure in Fe-based superconduc¬ 
tors. The most commonly discussed is the 1 S 4 = E 2 • C 4 
where a four-fold rotation C 4 is followed by a mirror re¬ 
flection about the x-y plane to account for the staggering 
positions of the Se/As ions above and below the Fe-plane. 
The spatial-orbital operations of S 4 is therefore 

S4 : x -A y,y -A -x,z -A- -z, 

d e -A —d e d e -Ad 1 d l -A -d t (A) 

^xz 7 ^yzi^yz 7 ^xzi^xy 7 ^xy i v*) 

d x 2_ y 2 -A -d x 2_ y 2,d z 2 -A d z 2. 


The rotation mixes the orbitals defined on the A and 
B sublattices and turns them into the corresponding 
even-odd combinations, d e J° = ± d^), that form 

the eigenstates of G s with eigenvalues: 4' tT (A:) = 
[d~(k), d+(k)] T where 


d~/+ — (d e /° d e !° d°^ e d°! e d°^ e ) 

— \ u xz<n u yzai u x 2_y2 a i ^xycn u z 2 a>' 


( 8 ) 


Since [H t (k),G s \ = 0 ,U block-diagonalizes H t (k) in the 
T a (k) basis, 


H k = UH t {k)U t 


P(k) 0 \ 

0 P{k + Q) J ’ 


( 9 ) 
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FIG. 2: (a) TB band of FeSe. Black solid lines are from P(k) 
and red dashed lines are from P(k + Q). (b) BZ of FeSe. The 
larger zone bounded by the black lines is the one-Fe BZ while 
the smaller zone with blue boundaries is the two-Fe BZ. Q 
corresponds to one of reciprocal lattice vectors in the two-Fe 
BZ. The dashed black lines with arrows indicate the trajectory 
over which the band dispersions are plotted in (a). 


where Q = ( 7 r, 7r) is the reciprocal lattice vector for the 
two-Fe unit cell. Eqs. © and © show that the lat¬ 
tice translation symmetry corresponds to their invariance 
under k — > k + Q, provided the exchange e P o is exe¬ 
cuted in orbital space [dill . This results in the important 
identification \d% a ,M') = \d° aa , M' + Q) = \d° a(T ,M) in 
the two-Fe zone. Diagonalizing P(K) gives rise to five 
band dispersions with —1 eigenvalue under glide opera¬ 
tion, and shifting them by k — > k + Q generates the other 
five with +1 eigenvalue under G s . These 10-bands are 
shown in Fig.2(a) for FeSe by the solid and the dashed 
lines respectively. Note that the band degeneracies at 
r are between the eigenstates within P(k) (odd under 
G s ) or within P{k + Q) (even under G s ), whereas the 
degeneracies at M are between the eigenstates of P(k ) 
and P(k + Q). 

On general grounds, the two-fold band degeneracies 
indicate that the space group of a Hamiltonian H must 
have at least one two-dimensional irreducible representa¬ 
tion. Specifically, there must exist at least two group gen¬ 
erators, representable by two matrices D\ and D 2 with 
[H, Di\ = [H , D 2 \ = 0, that are mutually noncommuting, 
[Z?i, D 2 ] ^ 0. For the noninteracting system described by 
the TB Hamiltonian H t , these two symmetries are man¬ 
ifestly the S 4 and E^. When interactions are included, 
the degeneracies will remain so long as a pair of D \ y2 
exist and the two-dimensional irreducible representation 
remains intact. For symmetry-breaking interactions that 
reduce the space group down to an abelian group with 
only one-dimensional representations, the band degen¬ 
eracies will be lifted as we will show below. 

Since the Hamiltonian Hk in Eq. © is manifestly di¬ 
agonal and the states in Eq. © are the eigenstates at T 
and M points [29| , the interplay between symmetry and 

band degeneracy can be studied by focusing on the dfj° 
and dy° orbitals at T and M points separately. Eq. © 
shows that the degeneracies at T are within the e-e and o- 
o pairs, i.e. \d e xz /d yz ,Y) and \d xz /d yz ,T ), whereas those 


at the M point, with the M' folded to M by the recip¬ 
rocal lattice operation discussed above, are between the 
e-o components | d xz /d yz , M) and \d%,/d° z , M) as shown 
in Fig. 1(c). These degeneracies clearly originate from 
the S 4 and the E^ symmetries defined in Eqs.© and © 
since, 

s 4 K„r> = -\d e yz ,r), s A \d e yz ,r) = \di z ,r), 
e d K„r) = \d e yz ,T), E d |d® z ,r) = \d e xz ,T), 

and similarly for the odd components at T; and 
SaK z ,M) = -I d° yz ,M), S 4 \d e yz ,M) = \d° z ,M). 

E dK„M) = | d° yz ,M), E d \d yz ,M) = \d° z ,M), 

and similarly when “e” and “o” are interchanged at M, 
which is equivalent to M' under the reciprocal lattice 
vector. More importantly, it is crucial that the two sym¬ 
metries do not commute, i.e. [S 4 , Ed] ^ 0 such that they 
form the two-dimensional irreducible representation of 
the space group necessary for the two-fold degeneracy. 

For a deeper understanding that will facilitate the 
classification of the possible degeneracy-lifting interac¬ 
tions, let’s recall the eigenbasis of the glide symmetry 
'F CT and Eq. ©. In the decoupled subspace of xz and yz 
orbitals,= (d e xza , d° xza , d e yza , d° za ) T . It is convenient 
to regard this four-spinor as the direct product of two 
two-spinors spanning the orbital space ( xz/yz ) and the 
sublattice space (e/o). In the rest of the discussion, the 
spin indices are suppressed for notational simplicity un¬ 
less otherwise noted. The TB Hamiltonian Hk in Eq. © 
can thus be written down explicitly up to constants at T 
and M points, 

U r = A r <7o ® t z , H m = Am<t z ® t z , ( 10 ) 


where eq and r* are Pauli matrices acting in the orbital 
and sublattice spaces, respectively; and Ar and A m are 
half the energy separations between the pair of degen¬ 
erate points at T and M. In this representation, the 
symmetry operators are 

(S 4 , E d , G s ) r = (-ivy ®r 0 , a x giro, a 0 ®T g ), (11) 
(S' 4 ,E ( j,G s )m = (~icr y <8 t x ,<t x ®t x ,<to®t z ). ( 12 ) 

They commute with H^/m correspondingly. 

Next, we show that the band degeneracies at T and 
M are protected by nonunitary symmetries in very much 
the same way that the time-reversal symmetry protects 
the Z 2 topological insulators Q. Note that the direct 
product of the Pauli matrices comprises the identity and 
15 generators. A given Hq = H^/m anticommutes with 8 
of the other 14 operators and commutes with the remain¬ 
ing 6. Remarkably, the 8 operators organize into two sets 
of 4 operators; each forms an independent Clifford alge¬ 
bra with Hq and is protected by the hidden antiunitary 
symmetry Xj or T 2 with Xj 2 2 = —1 respectively. We can 
therefore group the generators according to the following 



H 1 ,... 1 H 4 ,H 5 l ... 1 H 8 ,H 9 ,...,H 1A , (13) 


Tj : even even odd 

T 2 : even odd 


odd 

odd 


even 
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where [JTo, -ffg. 14] = 0, and 

{H 0 , = 0, T X H 0,1. 4TT 1 = 

(#0, #5,...,8} = 0, T2Ho^,...fiT 2 1 = #0,5,...,8, 

Tiffs,....uTT 1 = -ff 5 ,..., 14 , (14) 

T2-ffl ) ...,4 > 9,...14T 2 1 = —Tfl.4,9,...14- 

Thus, the two-fold degeneracy of the quantum states 
at F and M in f/o is protected by at least one nec¬ 
essary nonunitary symmetry T £ (Tj,!^), analogous 
to Kramers doublet (|a),T|a)) protected by the global 
nonunitary time-reversal symmetry. Specifically, we find 
that for Hr, 

T[ = ia y K (g) r 0 , Tj = ia y K 0 r z , (15) 

where K is the complex conjugation operator. They are 
effective time-reversal operators in the orbital angular 
momentum channel embedded symmetrically and anti¬ 
symmetrically into the sublattice space respectively. Op¬ 
erators that are even under Tf 2 are 

Tj r - even : H r , a 0 0 t x , a x 0 r y , a y 0 t v , a z 0 r y 
T 2 - even : H r , a 0 0 r y , a x 0 r x , a y 0 t x , a z 0 r x . 


Similarly for Hm at M point, 

T™ = a x 0 ir y K, T 2 m = ia y K 0 t x , (16) 

and 

T™ - even : H M , o x 0 r 0 , a y 0 r 0 , a z 0 r x , a, 0 t v , 
T 2 m - even : H M , a y 0 r z , <r x 0 r z , cr 0 0 r x , cr 0 0 t v . 

The identification of the antiunitary T-symmetries that 
protect the band degeneracies in Fe-based superconduc¬ 
tors is one of the main results of this work. The finding 
allows us to characterize the form of electronic order in¬ 
duced by the possible effective interactions with respect 
to their ability to lift the band degeneracy at T and M 
points. The experimental observation of the momentum 
space anisotropy, temperature dependence, and domain 
features associated with the degeneracy lifting can then 
be used to determine the important microscopic elec¬ 
tronic interactions responsible for nematicity in FeSe. 


J 


Interactions at T/M 

Real Space 

Order 

Ar 

Am 

Ql{{_^xz,edxz,e "b dx Z ^ 0 d X z,o') (^Z^)} 
92{_(dxz,edxz,e d) xz ^ 0 d xz ,o) 
9s{{dx Zje d xz ,e d) cz , 0 d X z,o') H - 

9^{d\. zi d X z,i ~ dy z i d yZ} i) 

9^{d\. zi d xz ,i+5 — d} yzi d yZl i+s) 

P 3 ( 1 ) v {d\. z ^dxz,i+S + d\ z ^dyz,i+s) 

Ferro-orbital 

s-wave bond nematic 

d-wave bond nematic 

0 

0 

0 

0 

0 

0 

94+(d xz ,edyz,e H - d xz , 0 dy Zi0 ) H - h.c. 

^9^+{.d^x Z , e dyz,e H - d xz o dy Zj0 ) ~b h.c. 
94—{d xz e dy Z ,e d XZj0 dy Z , 0 ) "b h.c. 
^ 4 -(^xz,e^j/z,e d^j. z ^ Q dy Z ,0) "b h.C. 

9^+^x Zt %dy Zl i ~\~ h.c. 

^9^-\-d\. zi d yz ,i + h.c. 
94:-{.d\ zi d yz ^5 + dy Zi d XZj i^s) + h.c. 
W4—(.dx Z jd yZt i+8 — dy Zi d XZ] i-\-s) + h.c. 

Orbital polarization 

Ferro Lz/Spin orbit coupling 
inter-orbital hopping 
Orbital current/Spin orbital flux 

0 

0 

0 

0 

0 

0 

0 

0 

9b+(A XZ ,edxz,o ~\~ dy Z e dy Z ^o) "b h.c. 
ig5+(d XZje d XZf0 ~b d\j Z e dy Z ,o) ~b h.c. 
9b—(d^xz,edxz,o dy Z e dy Zl o) ~b h.c. 

'i'95—(d xz , e d XZ)0 dy Z e dy Zj0 ) ~b h.c. 

p 5 _l_e z( ^ ri {d\. zi d xz ,i + dy Zi d yZj i) 
‘i95+{.d\ z ^d xz ,i+8 + dy Zi d yZj i-\-s) + h.c. 

95- e ^ ri {d} xzi d xz ^i — dy Zi d yZj i) 
^9b—{.d\ z ^d xz ,i-\-8 — d\ z idyz^i+s) + h.c. 

Charge/spin density wave 
circulating current (flux) 

AF orbital 

Ferro orbital current 

0 

0 

0 

0 

0 

0 

0 

0 

96+(.d X z,edyz,o "b d XZj0 dy Z ,e) ~b h.c. 

'i'96+(dt; Z ,edyz,o "b ^x2,o^J/2,e) "b h.C. 

96—(.d]cz,edyz,o d XZj0 dy Z , e ) “b h.c. 
i j 9Q—(.d xz ,edy Z ,o d XZ:i0 dy Z ie) ~\~ h.c. 

g's+e xQri d\ zi dy Z ,i + h.c. 
ig6+e lQri dl zi d yz ,i + h.c. 
ge—(—l) y {d xz i/ d yz ,i+5 — dy Z i d XZj i-\-s) + h.c. 
^96— ( 1 ) y {d>x Z 'id yZt i+5 ~b d\ Jz i d XZ} iJ t s) + h.c. 

AF orbital polarization 

AF L z l AF spin orbit coupling 
d-wave interorbital hopping 

d-wave orbital current 

0 

0 

0 

0 

0 

0 

0 

0 


TABLE I: The six types of effective interactions in the eigenbasis of glide symmetry (1st column), their real space representations 
(2nd column), physical meanings (3rd column), and whether they generate degeneracy splitting at T and M (4th and 5th 
column). 


B. Degeneracy lifting interactions 

On physical grounds, two degenerate states ipi and '02 
can be split by either level shift or quantum mixing (hy¬ 


bridization). Using the spinor notation, 0 = (ipi,ip 2 ) T , 
we can express the two types of interactions as gzi/j^aztp 
and g x il)'cr x xl>, respectively. The outcome of the a z in¬ 
teraction dependents on the sign of g z , which determines 
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the relative position of ip i and ip 2 upon splitting and 
hence corresponds to two possible domains. On the other 
hand, the a x interaction, being off diagonal, leads to a 
hybridized spectrum symmetric in ipi and ip 2 indepen¬ 
dent of the sign of g x , and hence produces no domain 
effect. Since ARPES experiments observed two-domain 
contributions for the splitting at M below T s and a split¬ 
ting with a single domain at T/Z above T s , they must 
originate from two distinct type of interactions, i.e. Am 
due to orbital shifting and Ap /z due to quantum mixing 
of the orbitals. 

It is therefore possible to write down all the possible 
forms of the electronic order or effective interactions be¬ 
tween the degenerate orbital states at T and M. Since the 
experimental observations are generally consistent with 
a unit cell containing two-Fe atoms, we will focus on the 
interactions that do not break the lattice translation sym¬ 
metry. We find six types of relevant interactions or elec¬ 
tronic orders among the d xz and d yz orbitals in the eigen- 
basis of the glide symmetry in Eq. flSJ). They are listed in 
the first column of Table I with the corresponding cou¬ 
pling constants < 71 ,..., ge- The coupling constants < 74 , 5,6 
are in general complex and are thus decomposed into 
real and imaginary parts by writing g m ± = g' m ± + ig'P n ±, 
m = 4, 5, 6 . In the second column of Table I, the spatial 
representations that produce the corresponding interac¬ 
tions at the T and M points are given. In cases where the 
interaction is nonlocal, we consider the nearest neighbors 
of site i indicated by * + 5 x ^ y along the x and y direc¬ 


tions and the lowest angular momentum representations. 
The physical meaning of the electronic interaction/order 
is given in the third column and whether it can lift the 
band degeneracy at T and M points is noted in the last 
two columns of Table I. 

In Table II, the matrix structures of the interactions 
are given explicitly in terms of the product of Pauli ma¬ 
trices in the orbital and sublattice subspaces. They span 
the complete set of the 15 generators. The columns in 
Table II show the properties of the corresponding inter¬ 
action under the point group and glide symmetry oper¬ 
ations given in Eqs m and m- Note that although 
these symmetry operators can have different forms at T 
and M in the latter equations due to the presence of the 
reciprocal lattice vector, the symmetry and symmetry¬ 
breaking patterns of the interactions < 7 ,; are the same at 
r and M, which is consistent with the fact that the space 
group operations are global. However, it is important to 
realize that these symmetry operations can have differ¬ 
ent commutation relations at T and M , thus resulting 
in different group properties such as whether there ex¬ 
ists a two-dimensional irreducible representation. Since 
the band degeneracy is protected by the newly identified 
antiunitary T-symmetry, it is straightforward to deter¬ 
mine if the interaction can lift the degeneracy by break¬ 
ing the Tf 2 and T-pfj at T and M points separately. In 
the case where the degeneracy remains, the protecting 
T-symmetry is displayed in the parenthesis. 



ffi 

CTz^~0 

92 

VzTz 

93 

croTz 

94+ 

OxTO 

// 

<?4+ 

CTyTO 

94- 

CfxTz 

// 

fl 4 - 

OyT z 

95+ 

O'OTx 

// 

<?5+ 

0~QTy 

95- 

O'zTx 

n 

95- 

(J Z Ty 

96+ 

&xT x 

II 

96+ 

&yT x 

96- 

OyTy 

II 

96- 

Gx'Ty 

54 

X 

X 

X 

X 

V 

X 

V 

V 

V 

X 

X 

X 

V 

V 

X 

Ed 

X 

X 

X 

V 

X 

V 

X 

V 

V 

X 

X 

V 

X 

X 

V 

G s 

V 

V 

V 

V 

V 

V 

V 

X 

X 

X 

X 

X 

X 

X 

X 

G s • 5 4 

X 

X 

X 

X 

V 

X 

V 

X 

X 

V 

V 

V 

X 

X 

V 

Gs ■ 

X 

X 

X 

V 

X 

V 

X 

X 

X 

V 

s/ 

X 

V 

V 

X 

Ar 

0 

0 

o(*) 

0 

0 

0 

0 

om r ) 

0(T 2 r ) 

0 (T 2 r ) 

o (tT) 

0(T 2 r ) 

0 (T 2 r ) 

0 (Tf) 

0 (T[) 

Am 

0 

o« 

0 

0 (T! M ) 

0 (T! M ) 

0(T 2 m ) 

0(T 2 m ) 

0(T 2 m ) 

0(T 2 m ) 

0 (T! M ) 

0 (T! M ) 

0 

0 

0 

0 


TABLE II: Symmetry properties of the effective interactions in the direct product space of orbital and sublattice. The laster 
two rows indicate whether the band degeneracy at T and M is lifted (0) or not (0). In the latter case, the corresponding 
protecting antiunitary T-symmetry is given in the parenthesis, or a * is given to indicate that the interaction is identically zero 
at the corresponding high symmetry point. 


1 . g 1,2,3 interactions 

Interactions <71.2,3 produce orbital shifts as can be seen 
from Table I and break both S+ and E d symmetries while 
preserving G s as shown in Table II. They are thus ne¬ 
matic interactions. Furthermore, since they break all the 
T-symmetries, the degeneracies at T and M are not pro¬ 
tected and the lifting of the degeneracy must come with 


domain effects. Thus, <71,2,3 are only suitable candidates 
for describing the observations at M points where two 
domains are observed by ARPES at low temperatures, 
(i) The real space expression of <71 is given in Table I, 

OfO = <7l ,dxz,ia — d yzilJ d yZt icr), (17) 

ia 

which coincides with the commonly discussed ferro- 
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orbital order (FO) parameter. It is isotropic in momen¬ 
tum space and leads to Ap = Am ^ 0 simultaneously by 
breaking all the T-symmetry, and is thus incompatible 
with the the experimental findings in FeSe. 

(ii) As shown in Table I, g 2 is a bond operator between 
the nearest neighbors. In momentum space, 

O s nb = 92 k(dl. z ka d XZy ka — dy z ka d yZ: ka), (18) 

ka 

where = cos k x + cos k y . This is clearly an extended 
s-wave nematic bond order parameter that breaks both 
S4 and E d symmetry. The s-wave form factor (7*,) van¬ 
ishes at M point, which makes it possible for a nonzero 
expectation of O s nb to generate a A r ^ 0 but A M = 0 . 
Although this splitting patten is consistent with ARPES 
at T > T s , O s nb cannot describe the high temperature 
isotropic phase since it would break the four-fold rotation 
symmetry of the Fermi surface at T with accompanying 
domain effects; both were not observed experimentally. 

(iii) In momentum space, the bond interaction 53 in 
Table I reads 

OdNB = 53 'y ( ldk{d} xz k a d XZt k(j + (19) 

ka 

where / 3 k = cos k x — cos k y . This corresponds precisely 
to the d- wave nematic bond interaction that describes 
the low-temperature nematic state. (OdNB} 7^ 0 leads 
to A m 7^ 0, but Ar = 0 since its form factor /?/. van¬ 
ishes at the zone center. Note that OdNB drives an in- 
phase d- wave bond order between the d xz and d yz or¬ 
bital, which should be contrasted to the out-of-phase 
symmetry-preserving d- wave bond between these orbitals 
already present in the hopping terms of the TB model 
j35l ]. In the next section, we will show how OdNB can 
generated by the intersite Coulomb interaction, resulting 
in a d- wave nematic state consistent with experimental 
observations at low temperatures. 

2 . 54,5,6 interactions 

The remaining 3 types of interactions, 54,5,6 in Table 
I, generate quantum mixing/hybrdization among the de¬ 
generate orbitals. Thus any resulting degeneracy split¬ 
ting would have only a single domain. Table II shows 
that although there are 4 interactions in each type with 
different space group symmetry properties, the degener¬ 
acy lifting pattern is the same within each type. 

(i) Interactions of the 54 type break either 64 or Ed 
symmetry while keeping the glide symmetry G s as seen 
in Table II. Thus, it is still possible for the remaining 
group to contain at least one two-dimensional irreducible 
representation. Indeed, one of the antiunitary symmetry 
in Ti 2 remains and protects the band degeneracy at M, 
whereas all T-symmetry is broken at T where the band 
degeneracy will be lifted as indicated in Table II. From 
Table I, the spin SU( 2 ) invariant representation of g r 4 , 


describes the orbital polarization due to an effective on¬ 
site crystal field correction, 

°OP = 54+ 'y 'X d]cz.igdyz,ia + h.C.). (20) 

ia 

Similarly, that of g’^_ generates an inter-orbital hopping 
or an extended s-wave orbital polarization, 

O s op = 54 - y ] r yk(d <! xz k(J dy Zy k<y + h.c.). ( 21 ) 

ka 

Although both Oop and O s op split the degeneracy at T, 
the fact that they both break S 4 and G s ■ S 4 symmetries 
makes them incompatible with the experimental obser¬ 
vation where the splitting at T in the high temperature 
phase maintains the four-fold rotation symmetry. 

Surprisingly, the imaginary components g'f± preserves 
S 4 symmetry despite of Ar 7 ^ 0 , which offers an intrigu¬ 
ing, and the only possible account of the observed proper¬ 
ties in the high temperature isotropic phase: degeneracy 
splitting only at T point, coexisting four-fold symmetric 
Fermi surfaces, and the absence of domain effects. The 
spin SU(2) invariant representation of g'f + is 

Ol* = *34+ YX^dyz,* — dy Zia d XZ} ia-), (22) 

ia 

which corresponds to an orbital angular momentum L z 
order that breaks the time-reversal symmetry. Remark¬ 
ably, there exists a time-reversal invariant but spin-SU(2) 
breaking representation 

Osoc = *54+ a ^xz,iadyz,ia ~ dy Z irJ d XZi i (T ), (23) 

ia 

which has an identical form as the spin-orbit interac¬ 
tion in the d xz / yz sector. Such an interaction can come 
from either the intrinsic SOC or be generated effectively 
0 . Since the experimentally observed Ar /z degener¬ 
acy splitting is nearly T-independent up to the highest 
measured temperature of 150K [H|,|T5}], we conclude that 
the latter originates from the intrinsic atomic SOC and 
estimate the strength of the SOC to be on the order of 
30meV. Noted that since the SOC involving all d-orbitals 
breaks the glide symmetry, the corresponding band cross¬ 
ings should be lifted except at M and A point in the BZ 
|29l l36j . We note that although not all such splittings 
have been detected at the present time, the hybridiza¬ 
tion between the d xy and d xz / yz orbitals near 1 /Z has 
indeed been observed by ARPES experiments [lj, [ll|. 

The analysis of the interaction g"_ in Table I can be 
made in the same spirit. The spin SU(2) invariant re¬ 
alization of this bond operator corresponds to interor¬ 
bital circulating current order (or the orbital flux phase) 
that breaks time-reversal symmetry. Similar to O soc in 
Eq. (1231) . there is a time reversal invariant realization 
that corresponds to spin-dependent orbital current order, 
where electrons with opposite spin-component traverse 
the lattice and accumulate opposite signs of the flux. 


(iii) Interactions of the type <75 and ge are hybridiza¬ 
tions between the even and odd orbitals and thus break 
the glide symmetry G s . The 55 interactions are diagonal 
in the orbitals. The real parts represent charge density 
wave or spin density wave order (g' 5+ ) at wave-vector 
Q and antiferro (AF) orbital density wave order (g' 5 _), 
while the imaginary parts are realizations of the orbital 
current (or flux) order (g"< ) and a hybrid of ferro-orbital 
and orbital current order (g as shown in Table I. Note 
that all individual interaction of the g§ type leaves one of 
the T-symmetry intact at T and M points as seen in Ta¬ 
ble II, and thus do not lift the degeneracy at T and M. 
However, it is interesting to note that the coexistence 
of pairs of interactions, e.g. g' 5 ± or g" ± , would remove 
all the T-symmetry protections and lead to simultaneous 
nonzero degeneracy splitting energies at T and M. 

(iv) The interactions of the gg type have the same or¬ 
bital content as those of the g 4 , but scatter between the 
even and odd components in the sublattice space. As 
shown in Table I, they generate AF orbital polarization 
(gfg + ), AF spin orbital coupling or AF orbital angular 
momentum L z (<? 6 +); interorbital hopping in the d -wave 
channel (g' 6 _), and d -wave orbital current or flux (g g_). 
The degeneracy at T is protected by one of the remaining 
T-symmetry Tf 2 , while it is lifted at M since all g g inter¬ 
actions are odd under T- 1 A f ! . It is important to note that 
although Ar = 0 and A m 7 ^ 0 , ge cannot account for 
the experimental observations at low temperatures since 
S 4 or G s ■ S 4 remains a symmetry which implies that 
the Fermi surfaces maintain four-fold symmetry without 
nematic distortions. 

Based on the systematic symmetry analysis in this sec¬ 
tion, we conclude that the distinct degeneracy lifting ob¬ 
served by ARPES experiments at T and M originate from 
two independent interactions: the spin-orbit interaction 
in Eq. E3l) responsible for Ar 7 ^ 0 already at high tem¬ 
peratures above T s and the d -wave bond nematic inter¬ 
action in Eq. (flT)l) that produces the low-temperature ne¬ 
matic phase with Am 7 ^ 0 . 


III. INTERSITE INTERACTION AND d-WAVE 
NEMATIC BOND ORDER 

In the rest of the paper, we focus on how interatomic V 
leads to strong band renormalization in FeSe and induces 
an emergent nematic order (Ocjnb) 7^ 0 in the ground 
state. To this end, we study the extended Hubbard model 

H = H t + Hu + H v , (24) 

where H t is the tight-binding (TB) model given in Eq. dTJ) 
and studied in the previous section. The intra-atomic 
interactions Hjj are given by the standard multi-orbital 


Hubbard model 

Hlj — U ^ + (U “./) ^ ] 17'ia'ni{3 (25) 

i,ot i,a.<(3 

J ^2 Si a • S ip + J ^2 

i,cz^f.3 

where U and U' are the on-site, intra- and inter-orbital 
on-site Coulomb repulsions and J is the Hund’s rule ex¬ 
change coupling with U = U' + 2 J. Note that when the 
Hamiltonian (12511 is used to describe the complete set 
of d-orbitals, J should be understood as an average of 
the exchange interactions of the t 2 g and the e g orbitals 
since the effects caused by the difference in the latter are 
usually small in a cubic system (37| Hlj • The extended 
interatomic Coulomb interaction is given by 

H v = V ^2 : nin 3 : 

<*.!> 

where the “normal-order” sign indicates that the direct 
Hartree term depending on the total density rn = Yl a rii a 
is subtracted, since that part of the interaction has been 
already included in the LDA. The same is true when 
treating Hjj in the HF theory [ 25 ], [26|. Thus, our 
treatment of interactions is in the same spirit as the 
LDA+U+V approach [H, lioj . The importance of the 
extended Coulomb interaction in Fe-based superconduc¬ 
tors has been emphasized previously J3(| with a focus on 
the properties associated with such p-d charge transfer 
metals. In the down-folded Fe-only model studied here, 
the interatomic interaction in Eq. (1261) is between the 
nearest-neighbor Fe atoms. 

A. Quantum fluctuations due to intersite V 

It is well known that the Fe-pinictides band structure 
is prone to a collinear SDW order that is also present 
in the multiorbital Hubbard model [2(| [4l[. Writing 

( c iaa c i0 °') = l[n a + <jm a cos (Qaf ■ ri)]S a pS a(T ', where 
Qaf = (7r, 0) and n a and m a are the density and spin 
density in orbital a, a nonzero m a is most easily obtained 
in the weak-coupling Hartree-Fock theory J4l| which is 
reliable when U is small. Since the electron correla¬ 
tion strength is comparable to the bandwidth of the d- 
electron complex, a complete description of the collinear 
SDW metal phase with realistic parameters would re¬ 
quire a strong coupling approach that takes into account 
the correlation effects nonperturbatively [261]. Here, we 
will carry out the Hartree-Fock theory for the extended 
Hubbard model in Eq. EH) using as effective parame¬ 
ters U = 1.4eV and J = 0.2eV. We will show that the 
quantum fluctuations induced by the inter-site V lead 
to a renormalized band structure where the Fermi level 
sits close to the vHS, which allows a weak-coupling ap¬ 
proach to capture the leading instability, the suppres¬ 
sion of collinear SDW, and the emergence of the d- wave 
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nematic bond order. The Hubbard interaction Hu in 
Eq. (E5l) is decoupled in terms of the self-consistent in¬ 
ternal fields in the charge and spin sectors [26|: A a = 
^(2U—5J)n—^(U—5J)n a and h a = 
where ( n,m ) = ^2 a (n a ,m a ). The self-consistently de¬ 
termined ground state indeed has Qaf-SDW order in 
the absence of the inter-site interaction V. To treat the 
quantum fluctuations beyond LDA, we decouple Hy in 
Eq. (l26l) in the hopping/bond channel, 

H v = -V Y (Xifd\ a d j0 + h.c. - \ X ?f\ 2 ), (27) 

(i,j),aP 

where xtf = and the spin index is suppressed 

for simplicity unless otherwise noted. In the presence of 
translation symmetry, the real nearest neighbor valence 
bond between any pair of orbitals (a, /3) can be decom¬ 
posed into the lattice harmonics of different angular mo¬ 
mentum in fc-space: 

Xij = + ZXpxVk + ‘ZXpvVk + 2 XdPk\, (28) 

fc 

where Xk and /3fc are the extended s-wave and d -wave 
form factors and = i sin k x ^ are the p- wave form 


factors. Comparing Eqs (l27l) and (l28l) to the tight- 
binding model in Eq. Q, it is clear that such quantum 
fluctuations amount to renormalizing the hopping inte¬ 
gral 

ttf = tg-VxZf (29) 

between the nearest neighbors. There are in fact two 
classes of contributions generated by the inter-site inter¬ 
action analogous to the situation in a general renormal¬ 
ization group analysis: (i) the corrections to the existing 
hopping parameters that maintain the lattice symmetry 
and thus lead to the renormalization of the band struc¬ 
ture, and (ii) the spontaneous generation of new and sym¬ 
metry breaking hopping channels. Correspondingly, we 
can write 

tt _ jrb.r. , rrs.b. 

tly — rly + tly . 

Since the full expressions for Hy r - and Hy h ■ are rather 
lengthy, we shall display explicitly the terms involving 
only the orbitals: 


Hy — Vy ' 2 [Xs7fc(dj z .k.d X z,k + d\ lz kdyZ'k) + XdPk(d\. z k d XZt k d yzk d yZ} k) + X P xVkd XZt kdxy,k 
k 

+ xlt'nldl^ k d X y,k + xt A lkdi y k d X y, k + h.c] + (terms involving e g orbitals)], (30) 

where \s = (Xs 1 * * * + X 22 )/2 and Xd = ( Xd 1 — Xa)I 2- The symmetry-breaking part is 

Hy = — Vo 2 [A a 7 k{d^ xz k d xz<k — d^ vz k d yz , k ) + A dPk(d^ xz k d xz ^ k + d) vz k d yz ^ k ) 

k 

+ (X^ /22 ?7fe + X^J^r\l)d\ z i yzk d xz i yzM + (xf 7 k + xfPk + xfxtk + X^V V * k)dl z ,k d yz,k 
+ (Xs 4/24 7fe + Xd 4/2A Pk)dl z/yZ}k d xy ,k + XpyVkdx Ztk d X y,k + X 2 p iVk d l z , k dxy,k (31) 

+ (xf x n X k + xtWk)dly tk d X y,k + xtPkdl y , k d X y,k + h-c. + (t erms involving e g orbitals)], 


where A s = (xV ~ xf)/ 2 and A d = (x^ 1 + xf )/2- The 
second term in Eq. (13111 is precisely the d -wave nematic 
interaction in Eq. Cl and the last term is its counterpart 
in the d xy channel. Note that we have denoted the inter¬ 
action strength differently as V and Vo in Eqs m and 
m- Although the bare values stemming from the micro¬ 
scopic Coulomb interaction are expected to the same, the 
effective interaction strengths Vo V due to the effects of 
orbital polarization, screening, and other orbital depen¬ 
dent contributions. To develop more physical insights, 
we will vary V and Vo independently around Vq/V = 1 
in the calculations. 


1. Renormalization of the band structure in bulk FeSe 

We first discuss how V renormalizes the band struc¬ 

ture by switching off Hy° ■ (i.e. setting Vo = 0). The 

evolution of the self-consistently determined low-energy 
band dispersions are shown in Fig. 3. With increas¬ 
ing V, the value of E v h — 350meV in the TB model 
at V = 0 shown in Fig. 3(a) is renormalized close to 
the experimental value of 25meV at V = 0.73eV in 
Fig. 3(b), and then to coincide with the Fermi level Ep 
at V = 0.763eV shown in Fig. 3(c). With increasing 

V, notice that the Dirac crossings of the d xy and d yz / xz 

bands also move up toward the Fermi level. Eventually, 





0 


4 = 

s 

© 

+ 

(g) © 

(h) 

© 

0 

© 


© 

_+ 



FIG. 3: Renormalized band structure (a-d) in unit of eV and 
corresponding FS (e-h) at V = 0 (a), 0.73eV (b), 0.763eV (c), 
and 0.85eV (d). Red dots indicate vHS at d xz / yz degeneracy 
point at M. FSs in (e) are 87.5% of their actual sizes. 


the vHS is pushed above the Fermi level which now cuts 
through the Dirac nodes, while the hole bands sink be¬ 
low Ep at V = 0.85eV as shown in Fig. 3(d), realizing 
an interesting state of a Dirac semimetal. Correspond¬ 
ingly, Figs.3(e)-(h) show the remarkable evolution of the 
four-fold symmetric Fermi surfaces: from the very large 
LDA hole and electron pockets in Fig. 3(e) to the ob¬ 
served small elliptical electron pockets with prominent 
quasi-ID character in Fig. 3(f), and then to four flower 
pedals at E v h = 0 in Fig. 3(g), and finally to four Dirac 
points in Fig. 3(h). 

Thus, the immediate consequence of the strong band 
renormalization in bulk FeSe is the dynamically gener¬ 
ated proximity of the vHS to Ep in the electron liquid 
phase. It is remarkable that for the range of V values 
where the vHS is within 25meV of Ef in the renormal¬ 
ized band structure, the self-consistent solution of the 
ground state is always an electron liquid state that is 
stable against collinear SDW and charge density wave 
order. In Fig. 4(a), the 2D density of states (DOS) of 
the renormalized band at E v n = 0 is plotted, showing 
its logarithmic divergence at Ep. Correspondingly, the 
calculated two-particle static susceptibility Xo(q) shown 
in Fig. 4(b) displays a sharp peak at q = 0 due to the 
vHS, which appears even stronger than the peak at ( 7 r, 0) 
which dominated the DOS in the absence of V. Thus, 
the electron liquid state is only unstable toward a q = 0 
nematic instability, analogous to the nematicity due to 
the vHS proposed for bilayer Sr-ruthenates [ 43 . Inter¬ 
estingly, we find that the vHS strengthens considerably 
when the hopping integral between the d xz / yz and d xy 
orbitals, described by the term —2it x 4 sin(k x / y )d+ z / yz d xy 
in the TB model, is decreased. The TB fit to LDA bands 
gives t l x 4 = 305meV. Reducing t x 4 by a factor of 3 leads 
to a significantly flatter d X yband near T /Z as observed in 
the ARPES experiments [l4( and to reduced curvatures 
of the d xz / yz band at the vHS point. Figs 4(a)-(b) show 
that the corresponding DOS and the q = 0 susceptibility 
are significantly enhanced, suggesting a much stronger 



FIG. 4: DOS (a) and static susceptibility (b) of V- 

renormalized band structure at E v h = 0 and corresponding 
degeneracy-splitting energy Am (c). Black lines in (a)-(c) for 
TB fit t x 4 = 305meV and blue line for f 44 = lOOmeV. (d) 
Am versus V around vHS for f 44 = lOOmeV; curves from 
bottom to top are for Vo/V = 0.6, 0.8, 0.9,1.0,1.1, and 1.2. 
The dashed black vertical line in (d) indicates where the vHS 
coincides with the Fermi level. 


nematic instability of the Pomeranchuk-type. 


2. d-wave Bond Nematic Order 

We next show that the q = 0 instability corresponds to 
the d -wave nematic bond order. To this end, we first set 
V = 0.763eV where E v h = 0 (Fig.3c) and switch on Vo in 
the symmetry breaking part of the Hamiltonian Hy b ■ in 
Ea. (l3Tl) . The self-consistent solutions show that among 
all the symmetry breaking terms, the leading instability 
occurs precisely in the d -wave nematic bond channel with 
nonzero A^ and 411 Eq. lRffl) . which has the largest 
form factor at M-point. The degeneracy splitting en¬ 
ergy Am between the d xz and d yz orbitals is calculated 
directly from the self-consistently determined eigenstate 
energies at M-point and plotted in Fig. 4(c) as a function 
of Vo. It serves as a quantitative measure of the degree 
of nematicity in the ground state. 

To explore the range of V over which the nematic¬ 
ity is controlled by the proximity of the vHS to Ep , we 
plot in Fig. 4(d) the Aju-map as a function of V for dif¬ 
ferent ratios of Vo/V. The distance of the vHS to Ep, 
i.e. -E v h, varies with V, and the dashed vertical line 
marks the location where the vHS sits at Ep, i.e. where 
E v h = 0. The most notable feature of Fig. 4(d) is the 
existence of plateaux nearly symmetrically distributed 
around E v h = 0 that grows in width and height with 
increasing Vo/V. For Vo/V ~ 1, the latter covers the ex¬ 
perimentally observed E v h — 25meV. Parallel results are 
obtained for £ 44 = 305meV where a similar Am requires 
a larger Vo/V. 
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In Fig. 5, we plot the band dispersion and the FS 
calculated self-consistently in the d- wave bond nematic 
ground state in a single-domain. In Figs. 5(a)-(b), the 
renormalized Ep of the electron liquid is at the vHS, i.e. 
E v h = 0. The d- wave nematic order removes the vHS 
by splitting the band-degeneracy around the Fermi level 
and produces the 4-fold symmetry breaking FS pockets 
as a result of the d- wave Pomeranchuk distortion. In 
Figs. 5(c)-(d), E v n = 25meV in the electron liquid phase, 
which corresponds to the observed value for the nomi¬ 
nally undoped FeSe. Figs. 5(e)-(f) are obtained for 1% 
electron doped Fei.oiSe which may be closer to the as- 
grown samples used experimentally 0, H. In the 
latter two cases, we find good agreements with the mea¬ 
sured dispersions for a single domain shown in Fig. 1(d) 
near the M point as well as the shape of the FS pockets 

IF 


It is important to point out that the d- wave nematic 
bond order and the ferro-orbital order are not orthogo¬ 
nal and mutually exclusive. They belong to the same 
space group as it is clear from Table II that the in¬ 
teractions g 1 , 2,3 have the same symmetry breaking pat¬ 
tern. As a result, both A d = J2k Pk[n xz {k) + n yz (k)} 
and A fo = ~ n yz(k)\ are none-zero in the 

d- wave bond nematic state, deriving from a momentum 
distribution function n a (k) that breaks the S 4 symme¬ 
try, unless additional particle-hole symmetry is present 
which is not the case in FeSe. The conventional Pomer¬ 
anchuk FS distortion is thus generalized to the case with 
orbital-lattice momentum coupling. This explains why 
local probes such as NMR have detected FO order in the 
low temperature nematic phase Hlzl- We stress that, 
in the present theory, the nematic transition is driven 
by the d- wave bond order that couples to the intersite 
Coulomb interaction and the FO order is induced para- 
sitically, since the energy lowering in going from the elec¬ 
tron liquid to the nematic state primarily comes from the 
d- wave nematic bond order. Indeed, since the FO order 
parameter Apo only couples directly to the local intra- 
atomic interactions but not to V, the degeneracy splitting 
energy Am = 8 VoAd — \(U— 5J)Apo in the HF theory is 
dominated by the contribution from A d , and even more 
so due to the large Hund’s rule coupling J in Fe-base 
superconductors |26l 1431 ] making U — 5 J much smaller. 
This explains the insensitivity of degeneracy splitting en¬ 
ergy A r at the zone center to the d- wave nematic bond 
order at low temperatures, which is estimated to have 
only about 5 meV variations from 120K down to 22K in 
the ARPES data [3], as shown in Fig. 1(d). Neverthe¬ 
less, the induced Apo does lead to the distortion of the 
hole FS near the zone center. For the values of U and 
J studied here, the hole FS pockets shown in Figs 5(d), 
5(f), and 5(h) near the F point indeed break the 4-fold 
rotational symmetry in a manner consistent with the ex¬ 
perimental observations [isj . 

In Figs 5(g)-(h), we show the results obtained when 
an atomic SOC term involving the entire 3d complex, 
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FIG. 5: Band dispersion in unit of eV (left panel) and FS 
(right panel) in d- wave nematic state, (a) and (b): V = 
0.763eV, E vH = 0, and Vo = 1.3eV. (c) and (d): V = 0.73eV, 
Ey h = 25meV, and Vo = 1.34eV. (e) and (f): V = 0.75 and 
Vo = 1.325eV at 1% electron doping, (g) and (h): With SOC 
Aaoc = 28meV, V = 0.763eV, Vo = 1.32eV at 1% electron 
doping. 


Hsoc = E A soc (a| L 1/3) (cr| S \a') df acr dip a ', is added 

iatfitjc t' 

to the Hamiltonian with A soc — 28meV. As predicted by 
the symmetry analysis in the last section, the SOC splits 
the band degeneracy at T without affecting the the bands 
at M. Moreover, it pushes one of the two hole bands 
at r below Ep, leaving a single two-fold symmetric hole 
pocket consistent with what was observed in ARPES and 
quantum oscillation experiments flU . ITFI ]. It is important 
to note that since the glide symmetry is broken by the 
SOC, the Dirac crossings located below the Fermi level 
between the d xy and d yz orbitals in a single domain are 
lifted by the SOC, as shown in Fig. 5(g). This small gap¬ 
ping of the Dirac points can serve as a landmark for the 
presence of a sizable SOC in FeSe, although the detec¬ 
tion of the d xy band by ARPES below the Fermi level 
has been notoriously difficult in bulk FeSe. In contrast, 
we expect that Dirac like crossings formed by the d xy and 
d yz bands originating from different domains will not be 
gapped. 

Finally, for all the cases studied with E v h within 
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25meV of the Fermi level, the d -wave nematic order dom¬ 
inates and the collinear SDW order is absent in the self- 
consistent solutions of the ground state. It is thus highly 
conceivable that magnetism and d-wave nematicity are 
competing caricatures in FeSe superconductors. In the 
present theory, the strong band renormalization, the sup¬ 
pression of the collinear magnetic order, and the emer¬ 
gence of the electronic nematic order have the common 
origin which is the inter-site Coulomb interaction. 


IV. SUMMARY AND DISCUSSIONS 


We have shown that the rise and the demise of symme¬ 
try protected degeneracies in the electronic band struc¬ 
ture can be used to probe the novel quantum states 
and the underlying interactions in correlated multiorbital 
electron materials. This work makes two advances in 
this direction with specific emphasis on Fe-based super¬ 
conductors. First, a systematic symmetry analysis re¬ 
vealed the “hidden” antiunitary T-symmetries that pro¬ 
tect the degeneracies at high-symmetry points in the BZ, 
and their connection to point-group and glide symmetry 
operations. This enabled the identification of the rele¬ 
vant electronic order/interaction that can break the T- 
symmetry and lift the band degeneracy. These results 
are applicable to all Fe-based superconductors. For bulk 
FeSe, the above analysis combined with recent experi¬ 
mental observation of the splitting of the band degen¬ 
eracy, including their momentum space anisotropy, tem¬ 
perature dependence, and domain effects lead uniquely 
to the conclusion that the splitting present already in 
the high temperature electron liquid phase at T is due to 
the atomic SOC, while the splitting at M is due to the 
d -wave nematic bond order that emerges only in the low 
temperature nematic phase. 

A microscopic theory is then developed to show that 
the unusually large band structure and FS renormal¬ 
ization, the absence of magnetism, and the emergence 
of the d-wave nematic order in FeSe can be explained 
by the important Fe-Fe interatomic Coulomb repulsion 
V. In addition to offering a natural description of the 
ARPES and quantum oscillation experiments, the the¬ 
ory suggests that the electronic nematicity is the driving 
force behind the tetragonal to orthorhombic structural 
transition. Interestingly, recent neutron scattering ex¬ 
periments in FeSe observed the spin fluctuations around 
( 7 r, 0) below the structural transition and the magnetic 
resonance at (7r, 0) below the superconducting transition 
[ 45 L 1 which can be explained by the presence of a ne¬ 

matic electronic structure ■ We also note that the 
same V-term has been argued recently to play an impor¬ 
tant role in stabilizing the s-wave pairing symmetry in 
Fe-based superconductors [44]. The overarching impor¬ 
tance of the extended Coulomb interaction V may origi¬ 
nate from the lack of the charge reservoir layers and the 
shorter Fe-Fe bond in bulk FeSe when compared to Fe- 
pnictides [23j . Interestingly, electronic nematicity with 



FIG. 6: (a) Schematics of the electronic structure near M 
showing the location of the vHS (F v h) before (a) and the de¬ 
generacy splitting energy Am after (b) the nematic transition, 
(c) Am as a function of E v h at fixed Vo/V = 1.8,1.7,1.6,1.4 
from the top line to the bottom line. 


similar phenomenology has been observed recently in 35- 
monolayer FeSe films with a larger A m = 80meV and 
higher T nem = 125K suggesting that further re¬ 

duced screening of extended Coulomb interaction in films 
can result in a stronger V and an enhanced nematic re¬ 
sponse. 

An important, falsifiable prediction of the present the¬ 
ory is the correlation between the emergence of the de¬ 
generacy lifting nematic state and the dynamical, inter¬ 
site Coulomb interaction V -induced proximity of the vHS 
to the Fermi level. It is thus desirable to seek direct 
experimental evidence for extended Coulomb interaction 
and further experimental tests by other techniques such 
as scanning probe and X-ray spectroscopy for the pres¬ 
ence of the vHS near Ep above and its removal below 
the nematic/structural transition. This state of affairs 
is summarized in Fig. 6 where A m is plotted versus the 
distance of the vHS to the Fermi level E v h at different 
values of Vo- It is in principle possible to tune the vHS 
by doping, pressure, or chemical substitution and study 
the corresponding changes in the nematic response such 
as the transition temperature and the band degeneracy 
splitting energy Am- While more experimental tests are 
clearly necessary, recent studies of chemically substituted 
bulk FeSe by S indeed find that A m increases with de¬ 
creasing E v h [50), qualitatively consistent with the pre¬ 
diction shown in Fig. 6(c). 

More importantly, since electron doping FeSe moves 
the Fermi level upward and away from the vHS at M 
point (see Fig. 6a), the present theory predicts that the 
nematic order will disappear while the vHS and band 
degeneracies survive at low temperatures when the ma¬ 
terial is subject to sufficient electron doping. Remark¬ 
ably, this has been observed recently by ARPES on bulk 
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FeSe whose surface layer is heavily electron-doped with 
Na [5lj]. The measurements show that the nematic state 
is absent and the vHS and the band degeneracy remain 
intact at 65meV below the Fermi level [511 ] . Surprisingly, 
a pairing gap near Ep develops at low temperatures that 
is consistent with the onset of a superconducting transi¬ 
tion at T c = 20K, much higher than the 9K transition in 
undoped FeSe, suggesting that the nematic state in bulk 
FeSe is a form of competing order of the superconducting 
state. Furthermore, in 50-monolayer FeSe films, a contin¬ 
uous reduction of the nematic order induced degeneracy 
splitting A m by surface electron doping with K has been 
observed, as well as an increase in the superconducting 
T c when nematicity is suppressed {52] . 

There is indeed an empirical correlation that the higher 
the electronic nematic transition temperature, the higher 
the optimal superconducting T c when electron doping re¬ 
moves the nematicity. Recent reports on 30-monolayer 
FeSe films, which are similar to the 35-monolayer films 
mentioned earlier with a nematic transition temperature 
around 125K {48{,|49|, shows a superconducting T c as high 
as 44K [52| under K surface doping. The strongest ne¬ 
matic phase with T nem ~ 180K is in fact an insulator 
observed in the nonsuperconducting N-phase of single 
and double-layer FeSe films grown on SrTiC >3 substrates 


[53l . [hi ] . Thermal annealing introduces significant elec¬ 
tron doping that removes the nematic insulating state 
in favor of the superconducting S-phase with the highest 
T c as much as 65K [5314551 ] . Our findings on the impor¬ 
tance of inter-site Coulomb interaction, the correlation 
induced proximity of the vHS near the Fermi level, and 
the new form of nonlocal, bond nematic orbital order 
with momentum space anisotropy provide considerable 
new microscopic insights into the intimate, competing 
relationship between nematicity and superconductivity, 
which may hold the key to understanding the pairing 
mechanism and to making T c even higher in these mate¬ 
rials. 
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